Statistical mechanics far from equilibrium: prediction and test for a sheared system 
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We report the complete statistical treatment of a system of particles interacting via Newtonian 
forces in continuous boundary-driven flow, far from equilibrium. By numerically time-stepping the 
force-balance equations of a model fluid we measure occupancies and transition rates in simulation. 
The high-shear-rate simulation data verify the invariant quantities predicted by our statistical theory, 
thus demonstrating that a class of non-equilibrium steady states of matter, namely sheared complex 
fluids, is amenable to statistical treatment from first principles. 
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Complex fluids relax slowly so their structure is rad- 
ically re-ordered by now, as in shear-aligning liquid- 
crystals jamming suspensions or liposome cre- 
ation 2| . Sheared fluids consist of particles following the 
same Newtonian equations of motion as at equilibrium, 
since no field is applied to drive them; only the boundary 
conditions differ. Nevertheless, they violate equilibrium 
statistical mechanics [f|, and only the distributions of 
entropy and work have been rigorously analysed in such 
cases (a, 0J Si El • I n processing and using complex fluids, 
a state of flux is the rule rather than the exception, e.g. 
molten plastic flowing into a mould, blood flowing within 
capillaries, or grease lubricating a rotating axle. Under 
continuous shear flow, these systems exhibit statistically 
steady states with intriguing similarities to equilibrium 
phase behaviour. For example, in "shear-banding" of 
worm- like micelles ff^, the fluid itself partitions the ap- 
plied shear into a region of low-viscosity oriented mate- 
rial at high strain-rate, coexisting with a slower, more 
viscous region. The parameters controlling this struc- 
tural phase transition are shear rate and concentration, 
in addition to temperature. Typically, simplified models 
with artificial dynamics (lfjj or near-equilibrium approx- 
imations [IH are employed in modelling these types of 
system, without knowledge of any fundamental princi- 
ples. Here we validate numerically a complete theory of 
the detailed statisitics of transition rates and occupancies 
in a realistic driven system with Newtonian interactions. 

We consider a macroscopic region of fluid, our sys- 
tem, embedded in a larger volume of the same fluid, 
which acts as a heat bath or reservoir, and exerts time- 
dependent random forces on the system's boundary that 
arc not predictable from a knowledge of the state of the 
system alone. An instantaneous microstate of the sys- 
tem is defined by the exact positions and momenta of 
all its constituent particles. The laws of motion govern- 
ing its dynamics can be fully summarized by a set {uj a b} 
of M transition rates between every possible pair of mi- 
crostates a and b that the system can adopt. Here, Lo a b is 
the probability per unit time that the system, currently 
in microstate a, will be found in microstate b an instant 
later (so uo a b = for transitions that would violate the 
laws of motion). In the presence of random impulses from 



the reservoir, the latter microstate is not uniquely deter- 
mined. Thus, the set {co a b} describes both the system's 
dynamics and the probability distribution of forces from 
the reservoir. 

Although forces from the reservoir are stochastic, their 
randomness is nonetheless governed by strict rules. The 
conditions of thermodynamic equilibrium constrain the 
transition rates to obey the principle of detailed balance 
(DB), which states that the ratio of forward to reverse 
transition rates between any pair of microstates must 
equal the Boltzmann factor of their energy difference, 
^ab/^ba = ex P(E a — -Eft) (with microstate energies Ei 
measured in units of the thermal energy ksT). In other 
words, the statistical properties of an equilibrium reser- 
voir impose M jl constraints on the M rates (one per 
pair). Equivalently, in a Langevin description of the equi- 
librium dynamics, the added noise must obey a Fluctu- 
ation Dissipation Theorem (FDT) [2l[ and is Gaussian, 
with strength determined by the thermodynamic tem- 
perature. In a sheared steady state, the system again re- 
ceives stochastic forces from the reservoir, but with some 
non-equilibrium distribution; indeed, those forces make 
the system flow. This generates a different set {u> a b} of 
transition rates, or equivalently a different noise distri- 
bution. 

In the absence of a rigorous theory of non-equilibrium 
statistical mechanics (notwithstan ding non-equilibrium 
generalizations of thermodynamics |22j), it has become 
common practice either to invent the non-equilibrium 
transition rates, or equivalently to assume that noise 
obeys the FDT or some other ad hoc criterion (such as 
colored, i.e. time-correlated noise). Although the con- 
cept of a non-equilibrium temperature is appealing and 
has been useful in interpreting non-equilibrium simula- 
tions (24j , statistical derivations typically make assumpa- 
tions about microscopic noise or rates [2{|. Since arbi- 
trary invention of the rates is forbidden at equilibrium 
by the constraints of DB, one might expect the same 
degree of constraint to arise also from the statistics of 
the sheared reservoir that influences the non- equilibrium 
rates. Indeed, a non-equilibrium counterpart to DB can 
be derived from either information-theoretic pol [27j or 
Gibbsian p8j arguments, which yield a one-to-one map- 
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ping between the set of rates {w°J} in the presence of 
an equilibrium reservoir, and those {w a b} for the same 
system (with the same Hamiltonian) bounded by the 
sheared steady-state reservoir. The existence of such a 
mapping means that the driven steady state has the same 
number of constraints as the equilibrium case, since one 
may first define an equilibrium system (respecting DB), 
then apply the mapping to derive the driven dynamics. 

For completeness, let us briefly discuss the ingredients 
required to derive the mapping [27l . |28| . The transitions 
are nontrivially correlated via the dynamics so, instead 
of microstate transitions, the basic objects for statistical 
analysis are phase-space paths (29| . Such a path describes 
a system's entire history of microstates; the position and 
momentum of every constituent particle at every instant 
during the very long duration of some steady-state ex- 
periment. A phase-space path T is therefore a complete 
description of the real physics exhibited by the system for 
a particular realisation of the noise, including all existing 
non-trivial spatial and temporal correlations. 

We conceptually construct an ensemble of weakly in- 
teracting systems that can (i) exchange energy (as in 
the equilibrium canonical ensemble) and (ii) arbitrar- 
ily distribute the total shear strain amongst the mem- 
ber systems (as occurs amongst a set of finite fluid ele- 
ments within a larger volume of fluid undergoing shear- 
banding). The crucial property of the phase-space paths 
of this ensemble is that they are uncorrclated with each 
other, except via those two quantities (energy and shear 
strain) that are exchanged over long range due to local 
conservation laws. Despite this lack of correlation, each 
path fully describes all the spatial and temporal corre- 
lations of its system. We consider all paths consistent 
with the laws of motion, but differing in their realisa- 
tions of the noise. For any such set of uncorrelated ob- 
jects, their exact probability distribution pr can be found 
[3l| by maximizing the entropy of that distribution (not 
to be confused with a thermodynamic entropy) subject 
to constraints from normalization and conservation laws. 
The derivation for the sheared case differs from that for 
equilibrium by only one additional constraint, fixing the 
ensemble-average of the total shear strain attained by 
the systems over the duration of the thought-experiment. 
The result [28[ is a relationship between the path distri- 
butions in the two ensembles, pr oc pp q expz/7r, where v 
is a Lagrange multiplier for the extra shear constraint, 
and 7r is the total strain for path T. By summing this 
relationship over all paths that contain a given transition 
between a pair of microstates, one obtains the exact map- 
ping, discussed above, between the sets of transition rates 
{w a b} and {w^}. (Supplementary material demonstrates 
how the same analysis can be used to obtain microstate 
occupancies, a non-equilibrium analogue of Boltzmann's 
law.) 

We have not stated the one-to-one mapping explicitly 
here, as it would require further notation to be intro- 
duced. However, it was recently noticed [13, that the 
mapping implies some remarkably simple relationships 



that, being exact, apply arbitrarily far from equilibrium. 
These relationships apply to a state-space with arbitrary 
connectivity between any set of microstates: 

1. The total exit rate from any given microstate dif- 
fers from its equilibrium value by a shear-rate- 
dependent constant that is the same for all mi- 
crostates, i.e. J2b ( u ab - Lull) = QM V a. 

2. The product of forward and reverse transition rates 
is the same in the equilibrium and sheared ensem- 
bles, i.e. UJab u ha = to c a l u) c b l V a, b. 

We have tested the above theory in a one-dimensional 
model system, a "fluid" of rotors (see Fig. QJi) , each in- 
teracting with its neighbours via torsional forces, and 
respecting Newton's laws of motion. The angular accel- 
eration of each rotor is proportional to its net unbalanced 
torque, which is the difference between the torques ap- 
plied by its two neighbours. The torque between neigh- 
bours has three contributions: conservative, dissipative 
and random. The conservative part is the gradient of 
the four-well potential U ( AO) shown in Fig. [TJa, which 
is a function of the angular difference between the ro- 
tors. Thus, the zero-temperature ground-state has all 
rotors parallel, but antiparallcl and perpendicular con- 
figurations arc also moderately favourable. The uncorre- 
lated random contribution to the torque (representing 
any microscopic degrees of freedom that are indepen- 
dent of shear strain, e.g. Brownian forces from solvent 
molecules) has a uniform distribution with width a and 
zero mean. The dissipative part is proportional to the 
difference in angular velocity between neighbours, with a 
constant of proportionality that can be rescaled to unity 
without loss of generality, leaving the noise strength a as 
the model's only parameter. Neighbours experience equal 
and opposite torques, so that angular momentum is ex- 
actly conserved in the model. The boundary conditions 
are periodic, and the equations of motion are numerically 
time-stepped. 

As described thus far, this is an equilibrium model. 
Once initial transients have died away, it exhibits Boltz- 
mann statistics in the occupancies of the potential 
U(A6). At sufficiently low noise strength, the relative an- 
gles between neighbours are mostly confined close to the 
local potential minima, with only occasional transitions 
between potential wells. The measured rates of those 
transitions respect DB. Note that DB is not imposed a 
priori; it emerges from the dynamics at equilibrium. 

We model a fluid under (angular) shear by twisting the 
model (see Fig.QJi). The twist is imposed via the periodic 
boundary condition, by introducing an offset in the an- 
gle measured between rotors either side of the boundary, 
and increasing that offset linearly in time, at a rate Lj, 
so that the twist rate per rotor is 7. As with ordinary 
periodic boundary conditions, this angular analogue of 
Lees-Edwards boundaries [H[ avoids edge effects, since 
it is a non-local condition on the topology of the space; 
the rotors cannot tell where the boundary is located, as 
they are only aware of relative angles. 
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FIG. 1: (a) The one-dimensional rotor model of length L. 
Each of the L rotors is characterized by its angle and an- 
gular velocity 9. (b) The potential of interaction between 
neighbours, U(A0) = — cos AO — cos4A6>, a symmetric peri- 
odic function of their angular difference AO, with four wells, 
labelled a, b, c, d. 



To apply our general theory, some region of the model 
must be defined as the system, while the large remainder 
is the reservoir, supplying unpredictable non-equilibrium 
forces to it. The system should be much larger than any 
correlation length, to ensure weak coupling to an uncor- 
related reservoir. Unfortunately, a large system implies 
a high-dimensional phase space (6 and 9 for each ro- 
tor), so that acquiring a statistically significant sample of 
all the transition frequencies becomes prohibitively time- 
consuming. We take two steps to reduce the phase space. 
First, we take the limit of small moment of inertia, so 
that momenta are no longer independent, and the phase 
space reduces to the set of inter-rotor angles Ad. This has 
the added advantage, in a one-dimensional force-chain, of 
reducing the correlation length to zero since, with van- 
ishing rate of change of angular momentum, the forces 
now balance globally. We are therefore able, secondly, to 
treat every inter-rotor gap (with its single characteristic 
variable A9) as a system, each surrounded by a non- 
equilibrium reservoir. Note that the general theory ap- 
plies to systems with non- vanishing correlation lengths, 
and we are treating a special case only for the sake of 
expediency. An impediment still exists for analysing this 
rotor model: the general theory applies to transitions be- 
tween microstates, whereas the measured transition rates 
are between potential wells a, b, c and d, that are finite 
in extent. If these four continuous sets of microstates 
are sufficiently analogous to true microstates, then the 
theory applies in this case. Subject to that qualification, 
we can use the model to test the theory's central assump- 
tion of ergodicity; i.e. the available phase-space paths are 
representatively sampled by the dynamics. 

We first test the predicted relationship between exit 
rates. The required quantity J2 a w ab implicitly depends 
on the unknown temperature of the compared equilib- 
rium system. However, we can eliminate that unknown 
by appealing to a symmetry of U (A9) . Since wells b and 
d are identical at equilibrium, they have equal total exit 
rates, + = + w^J. So relationship number 1 
predicts them also to have equal total exit rates in the 
driven case: LOba + ^bc = ^da + ^>dc for all imposed shear 
rates 7. This equality is not obvious, since the equilib- 
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FIG. 2: (a) Test of the prediction Uda + ^>dc = ^>ba + ^>bc V 7 
with noise strength a = 10. The left-hand ordinate mea- 
sures rates, using the same units as the abscissa, while the 
right-hand ordinate measures the dimensionless ratio, (b) 
Test of the predicted relationship ui a t uJta — ^>ad ^>da V 7, 
with a — 20. At this higher noise strength, higher shear 
rates are numerically accessible because of the greater number 
of observed backward transitions, (c) Test of the prediction 
Ucd^dc = LOcbUJbc V 7, for the same parameters as in (a), (d) 
Contours of the measured ratio (u>), a + u>bc) / (^da +^dc), pre- 
dicted to be unity across the whole parameter space of shear 
rate 7 and noise strength a. Data acquisition time limited 
simulations in the bottom left-hand corner. See supplemen- 
tary material for more data 



rium symmetry is broken in the driven case, where one 
of the potential wells is upstream of the other. Measure- 
ments of the four rates in question are plotted against 
shear rate in Fig. [5^ for a particular noise strength. The 
rates vary considerably with 7 and depart significantly 
from their equilibrium values. Nevertheless, the ratio of 
sums, as anticipated, remains very close to unity. 

Next we test the second predicted relationship, be- 
tween products of rates. Again, we exploit the sym- 
metries of the hypothetical equilibrium state to obtain 
a relationship between the measured rates in the actual 
driven system only. In equilibrium, symmetry of U(A9) 
(together with DB) implies = and = uj^. 
Substitution into proposed relationship 2 implies a con- 
straint on the measured rates in the driven system: 
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^ab^ba = u a du>da V 7. This prediction is tested in 
Fig. [2Jd: as before, while the individual rates vary sig- 
nificantly across the range of driving speeds, the predic- 
tion is obeyed to an excellent approximation. Similarly, 
the equilibrium symmetry about well c gives rise to a 
third non-equilibrium prediction, lu cc i LOdc = ^>cb ^>bc V 7, 
verified in Fig. [2t- 

To quantify the accuracy of the non-equilibrium the- 
ory across the model's whole parameter space, contours 
of the measured ratio {u>b a + u>bc) / i^da + UJ dc) arc plotted 
in Fig. The theory predicts a value of unity every- 
where. The small but significant discrepancies between 
theory and data may be due to the finite extents of the 
states a, 6, c, d, making them not true microstates. Al- 
though the non-equilibrium rotor model constitutes an 
imperfect test of the theory, the theory performs strik- 
ingly well here, and exceeds the predictive power of any 
approximate methods available to non-equilibrium statis- 
tical mechanics. Importantly, the discrepancies do not in- 
crease with 7, as would be the case for a near-equilibrium 
theory. The discrepancies at low a may alternatively 
be due to failure of the ergodic hypothesis. The theory 
shares this hypothesis with equilibrium statistical me- 
chanics, which also fails for non-ergodic systems such as 
glasses, but is successful for a wide range of applications. 



[1] L. A. Archer and R. G. Larson, J. Chem. Phys. 103, 3108 
(1995). 

[2] M. D. Haw, Phys. Rev. Lett. 92, 185506 (2004). 
[3] O. Diat and D. Roux, J. Phys. II (France) 3, 9 (1993). 
[4] O. Diat, D. Roux, and F. Nallet, J. Phys. II (France) 3, 
1427 (1993). 

[5] L. D. Landau and E. M. Lifshitz, Statistical Physics (But- 

terworth Heinemann, Oxford, 1980). 
[6] D. Collin, F. Ritort, C. Jarzynski, S. B. Smith, 
J. I. Tinoco, and C. Bustamante, Nature 437, 231 (2005). 
[7] J. Liphardt, S. Dumont, S. B. Smith, J. Ignacio Tinoco, 

and C. Bustamante, Science 296, 1832 (2002). 
[8] C. Jarzynski, Eur. Phys. J. B 64, 331 (2008). 
[9] C. Jarzynski, Phys. Rev. Lett. 78, 2690 (1997). 
[10] G. E. Crooks, Phys. Rev. E 60, 2721 (1999). 
[11] D. J. Evans, E. G. D. Cohen, and G. P. Morriss, Phys. 

Rev. Lett. 71, 3616 (1993). 
[12] G. Gallavotti and E. G. D. Cohen, Phys. Rev. Lett. 74, 
2694 (1995). 

[13] D. J. Evans and D. J. Searles, Adv. Phys. 51, 1529 
(2002). 

[14] S. M. Fielding, Soft Matter 3, 1262 (2007). 

[15] P. D. Olmsted, Rheol. Acta 47, 283 (2008). 

[16] A. Ajdari, Phys. Rev. E 58, 6294 (1998). 

[17] A. Santos, J. J. Brey, and J. W. Dufty, Phys. Rev. Lett. 

56, 1571 (1986). 
[18] T. Taniguchi and G. P. Morriss, Phys. Rev. E 70, 056124 



While the logistics of data acquisition has restricted 
our study of the rotor model to the zero-mass, zero- 
correlation-length limit, the theory should also apply to 
the more general case with momentum degrees of free- 
dom, thus encompassing phases with non-zero correlation 
lengths. However, even the testable model studied here 
exhibits highly non-trivial behaviour, excellently pre- 
dicted by the non-equilibrium statistical mechanical the- 
ory. Of course, no theory can describe all non-equilibrium 
steady states, since such states are much more diverse 
than the set of all equilibrium problems, for instance en- 
compassing molecular motors, convection cells, granular 
media and traffic flow. We have nevertheless presented a 
fundamental theory that governs the steady-state motion 
of any flowing system on which work is done by a weakly- 
coupled non-equilibrium reservoir that is ergodic (suffi- 
ciently mobile to explore its space of available states thor- 
oughly), and microscopically reversible (containing nor- 
mal particles with no sense of direction). Such systems 
include all sheared complex fluids, whose phenomenol- 
ogy is as important to future technologies as it is to our 
understanding of non-equilibrium physics. 
Acknowledgements: The work was funded by EPSRC 
grant GR/T24593/01. RMLE is supported by the Royal 
Society. 



(2004). 

[19] J. Casas-Vzquez, M. Criado-Sancho, and D. Jou, Euro- 

phys. Lett. 23, 469 (1993). 
[20] K. Kawasaki and J. D. Gunto, Phys. Rev. A 8, 2048 

(1973). 

[21] C. W. Gardiner, Handbook of Stochastic Methods 

(Springer, Berlin, 1985). 
[22] Y. Oono and M. Paniconi, Prog. Th. Phys. Suppl. 130, 

29 (1998). 

[23] S. I. Sasa and H. Tasaki, J. Stat. Phys. 125, 125 (2006). 
[24] T. K. Haxton and A. J. Liu, Phys. Rev. Lett. 99, 195701 
(2007). 

[25] L. F. Cugliandolo, J. Kurchan, and L. Peliti, Phys. Rev. 

E 55, 3898 (1997). 
[26] R. M. L. Evans, Phys. Rev. Lett. 92, 150601 (2004). 
[27] R. M. L. Evans, J. Phys. A 38, 293 (2005). 
[28] A. Simha, R. M. L. Evans, and A. Baule, Phys. Rev. E 

77, 31117 (2008). 
[29] G. E. Crooks, Phys. Rev. E 61, 2361 (2000). 
[30] T. S. Komatsu and N. Nakagawa, Phys. Rev. Lett. 100, 

30601 (2008). 

[31] F. Mandl, Statistical Physics (John Wiley & Sons, Chich- 
ester, 1988), 2nd ed. 

[32] A. Baule and R. M. L. Evans, Phys. Rev. Lett. 101, 
240601 (2008). 

[33] A. W. Lees and S. F. Edwards, J. Phys. C 5, 1921 (1972). 



